Single and Paired Point Defects in a 2D Wigner Crystal 
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Using the path-integral Monte Carlo method, we calculate the energy to form single and pair vacancies 
and interstitials in a two-dimensional Wigner crystal of electrons. We confirm that the lowest-lying 
energy defects of a 2D electron Wigner crystal are interstitials, with a creation energy roughly 2/3 that 
of a vacancy. The formation energy of the defects goes to zero near melting, suggesting that point 
defects might mediate the melting process. In addition, we find that the interaction between defects is 
strongly attractive, so that most defects will exist as bound pairs. 



At sufficiently low densities, the ground state of a 2D 
electron gas is expected to be a Wigner crystal || . Quan- 
tum Monte Carlo calculations || predict Wigner crys- 
tal formation for r s > 37 at zero temperature where 
r s = awl 'io is the dimensionless parameter character- 
izing the ratio of potential to kinetic energy. QMC cal- 
culations also predict that disorder stabilizes the solid 
phase relative to that of the liquid j|] thereby shifting 
the melting boundary to higher densities, r s rj 10. In 
light of recent experimental observations pj-ffypyli] |l7| ] 
of a new conducting phase in 2D in the density range 
10 < r s < 37, there is renewed interest in the proper- 
ties of a dilute electron gas close to the Wigner crystal 
melting boundary. Yoon, et. al. ]15| ] have observed a 
transition at r s w 40 which they attribute to the melting 
of a Wigner crystal. 

In this work, we calculate the energies of point-like 
defects in a clean 2D electron Wigner crystal. The mo- 
tivation for the study of defects is two- fold. First, local- 
ized defects are present in a finite concentration at any 
nonzero temperature (they have even been speculated to 
exist at zero temperature a super solid). Secondly, the 
melting process in 2D can be influenced or even deter- 
mined by defect formation |pjj|| . We investigate the en- 
ergy of two kinds of point defects: a vacancy at one of 
the lattice sites or an interstitial centered in one of the 
triangular unit cells. We define N^ei ( the defect index) 
as the number of electrons minus the number of lattice 
sites. Using Path-Integral Monte Carlo (PIMC), we com- 
pute the energy to introduce defects into an N-electron 
crystal with — 2 < Ad f < 2. We find that the lowest 
lying defect excitations are centered interstitials. Simi- 
lar results have been obtained by Jain and Nelson in 
hexagonal columnar crystals and by Cockayne and Elscr 
p0| in the 2D Wigner crystal. We then show that for 
TVdof = ±2, the defect creation energy increases as the 
defects are pulled apart. This indicates that the ground 
state for two defects is a bound pair. 

The system being simulated is composed of electrons 
confined to two dimensions (2D) and interacting through 
a repulsive 1 Jr potential and immersed in a positive back- 



ground with the same density. The Hamiltonian is 
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where energy is in units of the Rydberg ft 2 /2ma§, ao is 
the effective Bohr radius and lengths are in units of the 
Wigner-Scitz radius aw = (np)^ 1 ^ 2 and Vq a constant. 
To understand the role of defects in well-characterized 
electron crystals, we focus on systems with r s > 40. 

At large values of r s , the exchange contributions to the 
energy are small [p|pl[|. For sufficiently large r s , r s > 75, 
we consider distinguishable electrons and neglect anti- 
symmetry but include it for r s < 75. However, to keep 
the system stable, we forbid particle exchange by enforc- 
ing the condition |r, — Sj| < 1.1a where a is the near- 
est neighbor distance and Si is the i th lattice site. Such 
a "tether" is realistic because in a quantum crystal ex- 
changes are rare and the wave function is peaked around 
the lattice sites. We have verified the independence of 
our results by varying and removing this constraint. In 
these calculations we assumed a hexagonal lattice, shown 
to be the stable structure [2l| ]. 

To study defect formation in a Wigner crystal, we use 
the PIMC method |22]j . In this method, the density ma- 
trix 

p(Ro,R M ;j3) = J ... J dR l dR 2 ...dR M -ip{R ,Ri]T) 
p{Ri,R 2 ;t)...p(R M -i,Rm;t). (2) 

for the quantum system is evaluated by sampling paths: 
{R ,Ri...Rm-i,Rm}, and R k = {n )fc , ...r„ )fe } and r i)fc , 
a bead, is the position of ith electron in the fcth time 
slice and r = /3/M with (3 = l/k B T. The action 
— Iti(p(Rm-i, Rm'i t )) is evaluated by first splitting the 
potential into a long-range part and a short-ranged part. 
The long-ranged part (which is slowly varying) is handled 
in the primitive approximation, while the short ranged 
action is the exact pair action of two electrons. To evalu- 
ate the 2NM-dimensional integrals in Eq.(0), wc employ 
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the Metropolis bisection sampling technique p3| . Ewald 
pij sums are used to calculate the long ranged potential 
energy and pair action pS. 

At densities with r s < 75 we used restricted path 
integrals to account for Fermi statistics and stabi- 
lize the crystal against melting. Only paths which re- 
mained in the positive region of the Slater determinant 
were allowed. Such a restriction is exact if the nodal 
surfaces of the trial density matrix are correct. For the 
nodes we used a Slater determinant of Gaussian orbitals 
exp(— c(jCi— Sj) 2 ) with the lattice sites and c taken from 
[ p"T[ with a ferromagnetic spin arrangement. Calculation 
of the multiple exchange Hamiltonian necessary to deter- 
mine the ground state magnetic ordering is in progress 
but WKB calculations indicate that ferromagnetism is 
most likely @. 
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FIG. 1. a) Hexagonal 2D Wigner crystal with a single 6-co- 
ordinated vacancy and a 3-fold coordinated centered intersti- 
tial (CI) defect, b) Pair vacancy and CI defects. 

The simulations were performed at fixed temperature, 
area and number of particles. The energy to create N de f 
defects (see Eq. (^)) in a system of N particles and Ni = 
N lattice sites with a density of p = N/A is 

AE dc[ = [e(N + N dci ) ~ e(N)} {N + N de{ ), (3) 

where e(n) is the energy per electron for a system con- 
taining n electrons, a density p and Nt lattice sites. The 
densities and temperatures studied were r s =40, 50, 75, 



100, and 200 and T=5.0, 2.5 and 1.25xlO- 5 Ry. We used 
120 and 340 lattice sites. We performed separated PIMC 
calculations with different number of particles and lattice 
sites rather than more efficient procedures where parti- 
cles are inserted or removed [ fL3[ . Such differential pro- 
cedures are difficult because of the combination of large 
relaxation of the lattice, the very large zero point motion 
and the antisymmetry. Our calculations are very time- 
consuming as the system gets larger because one needs 
high accuracy to obtain the energy difference but the 
direct method allows better control over the systematic 
error. 

To illustrate the types of defects of interest, we depict 
in Fig. [l] vacancy and centered-interstitial (CI) defects 
in a 2D hexagonal lattice. The initial configuration of 
the vacancy and CI defects have 6-fold and 3-fold coor- 
dination, respectively. In PIMC all particles are dynam- 
ical variables, so that the surrounding lattice relaxes and 
may change the overall symmetry of the defect. However, 
the constraint makes such relaxation local. We find that 
the defect energy depends weakly on temperature. The 
following discussion concerns the defect energies for the 
temperature 1.25xlO _5 Ry. This is significantly less than 
the melting temperature and the Debye temperature. 

Fig. |^ demonstrates that the formation energy for an 
interstitial is consistently lower than creation energy for 
vacancies in 2D Wigner crystals at all densities. At r s — 
75 we find that AEqi ~ 0.65A£? vac . Simulations with 
340 sites give agreement with the results for the smaller 
lattice. Jain and Nelson |l^] and Cockayne and Elser |2C| ] 
have obtained a similar result. 

For a Coulomb system without fermion exchange, the 
energy of a defect can be expanded as 

E D = cir; 1 + c 3/2 r s 7 3/2 + c 2 r~ 2 ... (4) 

The first term is the static potential energy of the defect, 
the second the harmonic energy of the defect. Cock- 
ayne and Elser |2(J] have done exact calculation of c\ and 
c 3 / 2 - Shown in the lower panel of Fig. ^ is the anhar- 
monic contribution defined as the excess energy beyond 
the harmonic calculation. Anharmonic effects lower the 
defect energies to approximately half the harmonic value 
at r s = 50. Cockayne and Elser find that the defect en- 
ergy vanishes for CI defects at r s = 15 ± 1 and r s = 9 ± 1 
for vacancies. We obtain C2(vacancy) ~ —2.8 ± 0.2 and 
C2(CI) = —1.85 ± 0.15. From these values we estimate 
that Ed vanishes for interstitials at r s = 35 ± 2 and 
r s = 29 ± 2 for vacancies. The proximity of the vanishing 
of the interstitial creation energy to the melting density 
is highly suggestive that interstitial defects play a role 
in the melting process. It is possible that the assump- 
tion of the ferromagnetic spin arrangement is stabilizing 
the crystal with respect to interstitials for 35 < r s < 40. 
We have not considered the role of other defects such as 
dislocations. 
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FIG. 2. Formation energy (top figure) for vacancy (circle) 
and CI (square) defects as a function of l/r s at a tempera- 
ture of 1.25xlO _5 Ry for a system of 120 lattice sites (open 
and solid simbols are for Boltzman and Fermi statistics re- 
spectively). The triangles are results with 340 lattice sites. 
The lines (solid for vacancy and dashed for interstitial) are 
from harmonic calculations Jic|| . The anharmonic energy is 
shown on the bottom panel. 



We have also investigated the properties of pair defects 
as a function of the spacing between the defects. The 
geometry used for the study of the interaction energy 
between two defects is shown in Fig. [l]b. Shown in Fig. || 
are the pair binding energies for vacancies and CI defects: 



E bp = AE dci {N dcf = 2) - 2AE dei (N dci = 1) 



(5) 



as a function of the separation between the defects. As 
expected, for large separations the binding energy goes to 
zero. For all separations studied, the pair-binding energy 
is negative, indicating an attraction between defects. As 
is evident from Fig. ^, the interstitial pair-binding be- 
comes slightly positive at r s rj 50 indicating a possible 
weak-unbinding near the melting transition. No such be- 
haviour was observed for vacancies, however. In two di- 
mensions, any attraction is sufficient for the existence of 
a bound state of two defects. At r s = 100, the binding 
energy for interstitials pairs is roughly 50i^T, for vacancy 
pairs 125K and the classical melting temperature is only 
24K (assuming the effective mass of the electrons and 
the dielectric constant are unity). The symmetry of this 
paired state will be determined by the magnetic order 



of the lattice and the relative exchange frequencies of 
the defects. For a 3D Wigner crystal, Moulopoulos and 
Ashcroft [^7| have proposed that a paired Wigner state 
is a possible intermediary to melting. While our calcula- 
tions illustrate that interstitials are strongly bound, they 
do not indicate that a lattice of such paired defects is 
stable. 
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FIG. 3. Binding energy of pair defects at r s = 100 and 
T=1.25x 10 _5 Ry, computed for systems of 340 lattice sites. 
The geometry used to study such defects is shown in Fig. (lb). 
Circles are for vacancies and squares for interstitials. 

As seen in Fig. |^ the binding energy of the defect is 
much greater than the melting temperature, hence point 
defects are mainly interstitial pairs, which are expected 
to obey Bose statistics. 
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FIG. 4. Binding energy of pair defects separated by a sin- 
gle lattice site as a function of l/r s at T=1.25xl0 Ry, com- 
puted for systems of 120 lattice sites. Interstitials correspond 
to squares and vacancies to circles. 



Let us consider whether one could have a Bose conden- 
sation of defect pairs. The defect pairs are thermally ac- 
tivated so that their density is piD = Po exp(— E^d /ksT) 
with po = 1/tt and E^d the formation energy of a pair 
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defect. The Brezinski-Kosterlitz-Thouless (BKT) tran- 
sition occurs at a density p2i?A* /ksT » 0.28 (this is 
scaled from the helium data [^6|). A solution to these 
two equations exists provided that E^d < 0.42A*. Here 
A* = h 2 /(2m*) < 2/r 2 s since the pair defects are heav- 
ier than two free electrons. Hence a necessary con- 
dition for the existence of a super-solid phase is that 
E 2 d < 0.84/r,. This condition is shown in Fig. 0. Our 
results indicate that the supersolid transition could only 
occur in the Wigner crystal very near melting and is likely 
to be precluded by the unbinding of pairs. 
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FIG. 5. Phase diagram. The solid line is an estimate of 
the stable region of the Wigner crystal. The dotted line indi- 
cates the energy of the paired interstitial, the dashed line the 
condition for the BKT transition, E^d < 0.84/rf . 



As mentioned above, we find that in a 2D Wigner crys- 
tal most defects will be bound pairs of interstitials. The 
conductivity is likely to be dominated by defect trans- 
port. This could explain the extreme sensitivity of the 
conductivity in the insulating phase |(| to an in-plane 
magnetic field which would cause unbinding of singlet 
pairs. Further simulations are needed to decide whether 
pairing persists in the melted phase as has been proposed 
by one of us plj . 
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